#include <iostream>
#include <cmath>
#include <purg.h>


int main ( int argc , char **argv )
{
  // Allocate an array of 100,000 generators, with ids 0, 1, ... , 99999
  purg::parallel_gen  **gen_list = new purg::parallel_gen * [ 100000 ];
  for ( size_t i = 0 ; i != 100000 ; i++ )
  {
    // Each generator has the same seed with different id
    gen_list[i] = new purg::parallel_gen ( 0 , 0 , i );
  }

  // Estimate pi to 3 digits using Monte Carlo with numbers
  // from each generator.  This generates 1,000,000 doubles
  double result = 0;
  for ( size_t i = 0 ; i != 10 ; i++ )
  {
    for ( size_t j = 0 ; j != 100000 ; j++ )
    {
      // Get a random coordinate in [0,1) x [0,1)
      double x , y;
      gen_list[j]->get ( x , 0 , -1. );  
      gen_list[j]->get ( y , 0 , -1. );  
      double distance_from_origin = sqrt ( x*x + y*y );
      if ( distance_from_origin <= 1. )
      {
        result += 1;
      }
    }
  }

  // Compute pi
  double pi = 4. * result / 1000000.;
  std::cout << "Pi ~ " << pi << std::endl;
  std::cout << "Error ~ " << fabs(3.14159265368979323846 - pi) << std::endl;
}
